Midgap states and charge inhomogeneities in corrugated graphene 
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We study the changes induced by the effective gauge field due to ripples on the low energy 
electronic structure of graphene. We show that zero energy Landau levels will form, associated to 
the smooth deformation of the graphene layer, when the height corrugation, h, and the length of the 
ripple, /, are such that /{la) > 1, where a is the lattice constant. The existence of localized levels 
gives rise to a large compressibility at zero energy, and to the enhancement of instabilities arising 
from electron-electron interactions including electronic phase separation. The combined effect of 
the ripples and an external magnetic field breaks the valley symmetry of graphene leading to the 
possibility of valley selection. 
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INTRODUCTION. 

The discovery of atomic thick graphene layers [l], and 
the measurement of their novel electronic properties [2[ 
have led to a great deal of interest. The details of the 
structure of these layers is not completely elucidated, al- 
though there is clear evidence that free standing sheets 
are not fiatjl], and layers placed on an insulating sub- 
strate seem also not to be flat, following the corrugation 
of the susbtrate[4, 5]. The height and width of the ob- 
served ripples in free standing samples are consistent with 
models which take into account the tendency of carbon 
ions to form different types of bonds [6]. 

The intrinsic curvature of the graphene sheets 0, SH, 
10,|ll| , the modulation of hoppings by elastic strains [l2|, 
13|, [iJl, and hybridization between the tt and a bands 
induced by curvature [l^ enter into the effective Dirac 
equation which describes the low energy electron states 
as a deformation induced gauge field. In the absence 
of scattering between the two inequivalent valleys in the 
Brillouin Zone, these gauge fields change the phases of 
coherent electrons, and can play a role similar to that of 
an effective magnetic field H, [13, H , changing; 
the transport properties of the system in the presence of 
low (real) magnetic fields. 

Single layer graphene also shows inhomogeneities in 
the electronic distribution near half filling 0, These 
inhomogeneities may be related to the existence of 
ripples [111, [H. 

In the present work, we study the changes induced by 
the effective gauge field due to ripples on the low energy 
electronic structure of graphene. These changes become 
important when the effective magnetic length due to the 
ripples is comparable to the ripple size. Then, zero en- 



ergy Landau levels |2l| can exist in regions where the ef- 
fective field changes slowly. These midgap states are not 
suppressed by the off diagonal disorder associated to the 
ripples [22], and lead to a peak in the density of states at 
zero energy. While clean graphene is a semimetal, with 
a vanishing electronic compressibility at half filling, the 
existence of ripples lead to a large compressibility at zero 
energy, and to the enhancement of instabilities arising 
from electron-electron interactions. Note that zero gap 
states are expected when the flux of the effective mag- 
netic field is larger than one in a given region [23]. It 
has also been shown that a sufficiently strong random 
gauge field leads to a divergent density of states at zero 
energy 124. 



We present estimates of the strength of the effective 
field induced by ripples of the sizes observed experimen- 
tally, and discuss simple models which illustrates the for- 
mation of Landau levels as function of the deformation 
of the graphene layer. We also analyze the combined 
effects of ripples and a real magnetic field, showing that 
the equivalence of the two valleys in the Brillouin Zone of 
graphene is broken. We finally discuss possible instabili- 
ties which can arise due to the enhanced electronic com- 
pressibility, including the possibility of electronic phase 
separation. The Appendix describes simple models where 
the phenomena discussed in the paper can be studied an- 
alytically. 
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FIG. 1: (Color online). Left: sketch of the ripples considered 
in the text. Right: Modifications of the nearest neighbor 
hoppings induced by the ripple. See text for details. 



EFFECTIVE GAUGE FIELD IN SINGLE LAYER 
GRAPHENE. 

Qualitative estimates. 

Fluctuations of order At of the hopping parameter t 
around an hexagon of the honeycomb lattice give rise to 
an effective flux, in units of e/i/c, through the hexagon 
of order ~ At/t. If the hopping varies smoothly over 
a distance then At ^ dt{a/l), where St is the overall 
modulation of and a is a distance comparable to the 
lattice spacing. Then, the total flux in an region of area 
P is ^ ^ {dt/t){l/a). This number is also an estimate of 
the number of Landau levels with n = within that area. 
Finally, the magnetic length associated with the effective 
field is Ib ^ \/ {t/6t){la), and the separation between 
the n = and n = ±1 Landau levels is ^y {tSt){a/l). 
If the origin of the modulation St is solely due to the 
strains induced by a ripple of height h and length then 
dt ^ dlog{t)/dlog{dc-c)h'^ /{o^l)^ where dc-c ^ a is 
the length of the bond between nearest carbon atoms. In 
graphene, the parameter /3 = d \og{t) / d \og{dc -c) is ~ 
2. Hence, we find that the flux per ripple, in quantum 
units, is: 

For h ^ 1 — 5nm and / ~ 50nm, we find dt/t ~ Sv^/vy ~ 
10-2 - 10-\ and > 1. 



One dimensional ripples. 

We analyze numerically the emergence of midgap 
states as function of the modulations of the hoppings 
using a simple model with the geometry sketched in 
Fig. [1]. This translational symmetry along the y axis 
greatly simplifies the calculations. The changes in the 
electronic structure induced by the ripple are determined 
by the local value of the gradient of the tight binding 
hoppings. More complicated patterns of ripples can be 
decomposed in regions described by an average gradient 
of the hoppings. Hence, the model studied here should 
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FIG. 2: (Color online). Low energy states induced by a ripple 
as shown in Fig.[T]. The average hopping is = 3eV. The 
width of the ripple is 1200a = 168nm. The modulations of 
the hoppings are: Top left, 6t/t = 0, top right, 6t/t = 0.02, 
bottom left, 6t/t = 0.04, bottom right, 6t/t = 0.02, and a 
periodic electrostatic potential of amplitude 0.02eV. 

have similar features to other structures, provided that 
the magnitude of the gradients of the hoppings are com- 
parable. 

The translational invariance along the y axis implies 
that ky is a good quantum number. Hence, the hamil- 
tonian can be reduced to a set of effective one dimen- 
sional hamiltonians, one for each value of ky, which de- 
scribe the hoppings between the rows shown on the right 
side of Fig.[T]. In the absence of modulations, the abso- 
lute values of these hoppings are t and 2tcos{ky\^a/2), 
were a is the distance between carbon atoms. The low 
energy states are centered around kyV^a = 27r/3 and 
kyy/3a = 47r/3. The modulation of the hoppings leads to 
the replacement: 

t ^ t\\{x) 

2tcos(^) ^ ^t^(x) cos2 {(j)) + At^ix) sin^ {(j)) (2) 

where (j) = {^/3kya)/2. An external (real) magnetic field, 
described by the vector potential Ay{x) = Bx leads to 
the replacement of ky by ky + eAy/c. 

We have calculated the bands of the ribbons shown in 
Fig.[T] using the modulation: 

t^{x) = 6t sin (^^^ (3) 

Results are shown in Fig. [2], with / = 1200a (800 unit 
cells). Periodic boundary conditions are used. The plot 
corresponds to one of the valleys. The levels in the two 
valleys, in the absence of an external magnetic field, re- 
main degenerate. 

The calculations are consistent with the qualitative es- 
timates made above. The number of midgap states is 
proportional to the range of ky values where they are 
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FIG. 3: (Color online). As in Fig. [2], with a magnetic field 
B = 10 T. Left: K valley. Right: K' valley. 



defined. The higher Landau levels are less well defined, 
in agreement with the the fact that the description of a 
hopping modulation as an effective gauge field becomes 
exact only at the Dirac energy. The results are similar to 
those found in the analysis of the Quantum Hall Effect 
in thick nanotubes in a perpendicular field [25[, where the 
(real) magnetic field is modulated and has zero total flux. 
There is a region in momentum space where the gaps be- 
tween subbands have a minimum, corresponding to the 
region in real space where the effective field vanishes. 

The effect of a uniform external magnetic field, which 
has the same sign for the two valleys, is shown in Fig. [3]. 
The combination of the breaking of inversion symme- 
try, induced by the ripple, and time reversal symme- 
try, induced by the magnetic field, leads to the inequiva- 
lence between the two valleys. The total effective field is 
greater in one than in the other, as shown by the larger 



J 



region occupied by midgap states, and the fact that the 
higher bands show less dispersion. The results presented 
in Fig. [3] show that pseudomagnetic fields created by 
the ripples broaden all Landau levels except the zero- 
energy one. Interestingly, measurements of the activa- 
tion gaps in the Quantum Hall Effect regime in graphene 
demonstrate that zero-energy Landau level is much nar- 
rower than the other ones [26(. The results obtained us- 
ing this model are in agreement with calculations using 
the full valence band of graphene and the Local Density 
Functional Approximation for one dimensional ripples of 
smaller sizes 1271. 



Two dimensional ripples. 

We consider now a single two dimensional ripple with 
axial symmetry. The elastic strains in circular coordi- 
nates are: 



drUr + (drh) 

dgUr uq drhdeh 



uee 



drUe 



fy* fy* (y* 

The Dirac equation in radial coordinates is: 



(4) 



ie'%F ( + — ) ^A(r, 0) + tf3e-^'^ {urr - uee - iuro) ^A(r, 0) = e^sir, 0) 



(5) 



r 



where /3 = 9 log (t) /Slog (a), t is the hopping between 
nearest neighbor orbitals, and a is the distance between 
carbon atoms. If the gauge field induced by the ripple 
has circular symmetry, we obtain: 



t/3 {Ur 



Uee + iure) 



fir) 



(6) 



A sketch of the ripples studied here is shown in Fig. [4], 
and also the associated effective magnetic field. The pa- 
rameters used in the figure, / 600A and h = 30A give 
a fiux of effective magnetic field of order unity per ripple. 

We compute numerically the electronic levels of ripples 
with axial symmetry, embedded in hexagons with peri- 
odic boundary conditions. The hexagons are labeled by 
an integer such that the total number of sites that 
they contain is 6{N * 1)^. The side of the hexagon is 
L = y/3{N -\- l)a, where a = 1.4A is the distance be- 
tween carbon atoms. The strains are simulated by an in 



plane deformation which modulates the hoppings: 



with: 



and: 



tfiVi = to -\- toP- 



Sfi 



(7) 



(8) 



(9) 



The resulting strain tensor is the same, in the continuum 
limit, to that is induced by a gaussian modulation of 
the height of the graphene sheet, hif) = hoe~''^^^^ , with 
A = p{ho/l)\ 

In the absence of the ripple, the levels have usually six- 
or twelvefold degeneracy, and are separated by gaps of or- 
der vy/L. There is a fourfold degenerate state at £^ = 0, 
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FIG. 4: (Color online). Right: Sketch of the ripple with 
an axially symmetric profile, similar to those analyzed in the 
text. The length is / = 600 A, and the height, is /i = 30 A. Left: 
Effective magnetic field generated by the ripple, see text for 
details. 
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FIG. 5: (Color online). Energy levels (in units of t) for ripples 
with circular symmetry embedded in an hexagon with 6 x 
52^ = 16224 atoms and periodic boundary conditions, and 
different radii / (in units of a). The length of each line is 
proportional to the degeneracy of the level. Top: A = 0.2. 
Bottom: A = 0.4 (see text for details). 

as the number of unit cells in the hexagon is a multiple of 
three. The electronic spectra at low energies for hexagons 
with = 76, (the number of sites is 34656) with embed- 
ded ripples of different sizes are shown in Fig. [5]. The 
side of these hexagons is L ^ 131a ~ 18nm. The ripple 
is defined as in eq.([9|), with A = 0.2 and A = 0.4, and 
the radius of the ripple is / = 50a, 60a and 70a. The 
modified hoppings range from 0.9t to 1.2t for A = 0.2, 
and between O.St and lAt for A = 0.4. 

In the presence of ripples, the degeneracy of the lev- 
els is reduced, and the number of low energy modes in- 
creases with the size of the ripple, in agreement with the 
analytical estimates made above. Note that there are 
still significant gaps between these states, indicating that 
they will lead to peaks in the density of states, and not 
to resonances within a broad continuum. The periodic 
boundary conditions imply that the low energy levels of 
a given ripple are hybridized with those in neighboring 
ones, leading to a shift from E = 0. The dependence of 



FIG. 6: (Color online). Low energy levels of the ripple with 
A = 0.4 and / - 70a shown in Fig. [5]. 
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FIG. 7: (Color online). Top: Log- log plot of the scaling, as 
function of the size of the hexagon, of the lowest energy state 
away from ^ = of a system with a ripple of radius / = 50a, 
and amplitude A — 0.4. Bottom: Scaling of the fifth level 
away from E = 0, and for the same parameters. The red 
straight line has a slope of -1. 

the level spacing with 6t and / agrees with the estimates 
in Section Ila. 

The low energy spectrum for a ripple of radius / = 70a 
and amplitude A = 0.4 is shown in Fig. [6]. The level 
spacing near E = 0^ AE ^ 0.02t implies an effective 
magnetic length, Ib = v-p/AE = {3ta/2)/ AE ^ 75a ~ I. 
This length is larger, although the order of magnitude 
is comparable, than the analytical estimate obtained in 
Section Ila, h - ^{t/5t){la) - v'70/0.3 - 15a. Given 
the uncertainities associated to the estimate in Section 
Ila, and the fact that we have not taken into account 
that the effective field induced by the ripple is split into 
six peaks of smaller size than the ripple itself, we consider 
the agreement reasonable. 

The scaling of the lowest level close to £^ = as func- 
tion of system size, L, for a fixed radius of the ripple, 
/ = 50a, is shown in Fig. [7]. The energy scales approxi- 
mately as for / < L, and it shows a weaker depen- 
dence on L for / <C For comparison, the scaling of 
a level further away from E = {) shows the scaling 
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FIG. 8: (Color online). Top: Density profile of the wavefunc- 
tion of the state closest to ^ = 0, of a ripple with radius 
/ = 50a. The energy of this state is shown in Fig. [7]. Bottom: 
density for the fifth level away from ^ = 0. 




as function of system size expected from a delocalized 
level described by the Dirac equation. The radial distri- 
bution of the density associated to these wavefunctions is 
shown in Fig. [8]. In agreement with the scaling behavior 
of the energy, one of the states is localized within the 
ripple, while the other is extended towards the edges of 
the hexagon. 

All the lattices studied, with and without the modu- 
lation of the hoppings, show four zero energy states. In 
the clean system, the existence of these states is deter- 
mined by the valley and sublattice degeneracy, and the 
density associated to them is uniform throughout the lat- 
tice. The spatial extent of these states changes qualita- 
tively in the presence of a ripple, as shown in Fig. [9]. The 
wavefunctions become localized within the ripple, leading 
to a peak in the local density of states at £^ = 0, even if 
these states are delocalized. The wavefunction shows the 
hexagonal symmetry of the undelying lattice, and most of 
the charge is localized in six regions of the ripple. At the 
neutrality point there are only two (four including spin) 
electrons available in these states, opening the possibility 
of charge fractionalization. 

The results presented above have been obtained for 
large modulations of the hoppings, or, alternatively, 
strong height corrugations. This limitation is imposed 
by the dimension of the hamiltonian which can be diag- 
onalized. 

The calculations show, however, that the low energy 
part of the spectra have the features expected from the 
continuum Dirac equation, and the wavefunctions asso- 
ciated to these states have a smooth envelope on scales 
comparable to the lattice constant. Hence, the method 
used can be considered a finite element technique which 
approximates well the Dirac equation in the presence of a 




FIG. 9: (Color online). Top: Density as function of distance 
to the center of the wavefunction of the state at ^ = 0, of a 
ripple with radius / = 50a. Center: Density for sublattice A. 
Bottom: Density for sublattice B. The color code is such that 
density decreases from red to green, blue, and light blue. 



ripple. The property of the lattice model which needs to 
be kept invariant is v-p = (3ta)/2. The continuum equa- 
tions are invariant under the scaling / A/, St St/X^ 
and En En/X. The scaling of St implies that the cor- 
rugation of the ripple scales as (h/l)'^ {h/l)'^/X. As 
a result, we can extrapolate the results analyzed here to 
ripples of larger size and weaker corrugations. The result- 
ing effects are in agreement with the qualitative estimates 
made in Section Ila. 



INTERACTION EFFECTS. 

The analysis in the preceding section shows that rip- 
pled graphene samples, in the absence of interactions, 
have a peak in the density of states at the Dirac energy. 
The width of this peak decreases as e-(^^/^)(Vo^)^ and the 
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FIG. 10: (Color online). Sketch of the free energy of an or- 
dered phase induced by interactions. The ordered phase is 
assumed to have a gap, and the derivative of the free energy 
at n = is discontinuous. On the other hand, the disor- 
dered phase is assumed to have a high electronic compress- 
ibility. The Maxwell construction indicated by the straight 
lines shows the region where electronic phase separation oc- 
curs. For comparison, an ordinary second order phase transi- 
tion between phases with similar electronic compressibilities 
is also shown. See text for details. 



fraction of states, in an area of size that it includes 
is proportional to the same dimensionless combination, 
/ ~ {5t/t)(l/ a). The resulting diverging electronic com- 
pressibility implies that interactions will induce insta- 
bilities, in the same way as in graphene in a magnetic 



field near half filling |28|, [H, lH, [H, [H, S, [35[. As 
shown in Fig. [2], while the midgap states are well de- 
fined, the higher bands, for small values of 8t/t^ tend 
quickly towards the clean graphene limit. If the inter- 
actions are weak, this implies that their effect is limited 
to the states within a narrow range of energies near the 
neutrality point. Alternatively, it can be argued that 
interactions can be studied in a restricted hamiltonian 
which includes only the midgap states. Hence, a rea- 
sonable upper bound to the gap opened by many body 
effects is A - {t6t) {a / 1) . 

A plausible ordered phase which can emerge at low 
temperatures is a ferromagnetic phase [30., The satu- 
ration magnetization will be small, m ~ {6t/t)/{al). The 
Curie temperature needs not be small, as transverse spin 
waves in a system like this do not induce a large sup- 
pression of Tc [36[. This phase may help to explain the 
observation of ferromagnetism in graphite samples [37]. 

Phase transitions in systems where the electronic com- 
pressibility is large are often of first order, and the system 
has a tendency towards electronic phase separation in the 
ordered phase The expected dependence of the 

free energy of the ordered and disordered phases as func- 
tion of electronic concentration is sketched in Fig. [10]. 
The region of phase separation which can be obtained 
from a Maxwell interpolation using the free energies of 
the different phases will be replaced by a disordered phase 
with charge puddles, as the long range Coulomb inter- 
action suppresses macroscopic phase separation. The re- 
sulting situation resembles that in magnetic systems with 
striped phases or domains [40i, i41i]. 



We can make a qualitative estimate of the sizes of the 
charge puddles by comparing the energy gain per unit 
area in the ordered phase and the electrostatic cost of 
forming the puddle. For a puddle of size Lpud^ the energy 
in the ordered phase decreases by an amount AnL^^^, 
where n is the density of electrons contribute to the or- 
dered phase. The electrostatic energy is e^n^L^^^. Then, 
the typical puddle size is Lpud ^ A/(e^n). Putting to- 
gether the estimates made above, we obtain: 



-'pud 



As in graphene 




(10) 



0(1) for rea- 
sonable ripple parameters, the size of the puddle will 
be comparable with that of the ripple. A sketch of 
the gap opened by an periodic electrostatic potential, 
vix^y) = 'Uo sin(27rx//), with vq = 0.02eV is shown in 
the bottom right graph in Fig [2]. The gap is of the same 
order of magnitude as the amplitude of the potential. 



CONCLUSIONS. 

We present analysis of the changes in the electronic 
structure of graphene due to modulations in the hop- 
pings induced by ripples and other sources of elastic 
strains. The changes in the electronic structure are deter- 
mined by the dimensionless parameter ^ ~ {6t/t){l/a) ~ 
{Svj;^/vY){l/a) ~ {[3h^)/(la)^ where 5t/t = 6vy/vy is the 
modulation of the hopping parameter, or the Fermi ve- 
locity, vp, ^ is the height of the rippls, / is the size of 
the ripple, and a is the lattice constant. The parameter 
^ gives the fiux of effective magnetic field threading an 
area of the size of the ripple. 

We find that reasonable values for the size and height 
of a ripple lead to the formation of midgap states. These 
midgap states are similar to the Landau levels at the 
Dirac energy induced by a magnetic field. The combina- 
tion of hopping modulations and a magnetic field breaks 
the symmetry between the two graphene valleys, lead- 
ing to the possibility of valley selection [42*], as electrons 
from each valleys will scatter differently from extended 
defects jial. 

Midgap states induce a large electronic compressibility 
when the Fermi energy is near the Dirac point. Inter- 
action effects will lead to instabilities towards ordered 
phases, and electronic phase separation, with typical 
puddle size not too different from that of the ripple. 
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APPENDIX. ANALYTICAL MODELS OF LOW 
ENERGY STATES IN RIPPLED GRAPHENE. 

Straight ripple. 

We analyze simple models of elastic deformations and 
the midgap states that they may induce. 

We study first a straight boundary between a stretched 
and a relaxed region in graphene, as schematically shown 
in Fig.[TT]. We describe the change in the stretched re- 
gion by a change in the nearest neighbor hopping, t, along 
the horizontal direction, t -\- St. 

The Dirac hamiltonian in the stretched region, is: 



n = 



(11) 



v^i—ikx T ky) + St 

VF{kx T ky -\- St) 

so that the perturbation induces a gauge field in the y 
direction. Ay. At the boundary, we have dxAy ^ 0. 

The system has translational symmetry along the y 
axis, so that ky is conserved. The Dirac equation, for a 
given valley reads: 

-vyOxiI^b^x) ± {vpky T St{x))iljB{x) = eijjA{x) (12) 

These set of equations give the one dimensional Dirac 
equation with a gap A{x) = v-pky -\- St{x). If the gap 
changes sign across the interface they have localized so- 
lutions at e = 0. This condition implies that either 
v^ky > and v^ky^St < or v^ky < and v^ky^St > 0. 
Irrespective of the sign of St^ there is a range of values of 
ko,'. 



vf 



(13) 



The electronic structure of the system is schematically 
shown in Fig. [12]. The distortion shifts the Dirac cone in 
the stretched region. A band of localized states joining 
the two Dirac cones is induced. The number of midgap 
states per unit length of the ripple is: 

S\St\ 



niD 



3ta 



(14) 



where t is the nearest neighbor hopping, and a is the 
C — C distance. 




FIG. 11: (Color online). Sketch of a boundary between a 
stretched and a relaxed region of graphene. 




FIG. 12: (Color online). Electronic spectrum as function of 
ky of the interface shown in Fig. [11]. 



Circular ripple. 

The zero energy wavefunction in a given valley, in the 
presence of a circular ripple, obeys the Dirac equation, 
eq.©: 



ie^%F [Or + ^ J ^A(r, 0) + f{r)e-'''^A{r, 0) = 

(15) 

and related equations for ^ b (^, 0) and for the other val- 
ley. This equation can be integrated analytically, and we 
obtain: 

^(r, 0) = f{re'^)e'^ dr'r'^f{r') ^^^^ 

where f{z) is an analytic function. This wavefunction 
is not normalizable, irrespective of the choice of f{z). 
Choosing f{re^^) = A, where A is a constant, the elec- 
tronic density is: 

p{r,0)=A^e ^F^^ ^0 ^^"^ (17) 

Assumming, for instance, that lim^^oo /(^) ~ 
O ^e"^*^^/^^^^ ,and that /(r) does not diverge as r ^ 0, 
we find: 



lim p(r, (9) A 



lim p(r, 0) 



2 



(18) 
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The density has maxima at for a radial coordinate com- 
parable to the radius of the ripple, r ~ /, where it depends 
exponentially on a quantity of order f{l)l/vY ^ I3h^ / {la). 
As a function of 0^ p{r ^1^0) has three maxima and three 
minima. The positions of the maxima and the minima are 
interchanged when analyzing the other valley. The re- 
sulting pattern is in good agreement with the numerical 
results shown in Fig. [9]. 
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